Today

First hour

  • Hello, I’m healthiar!

  • healthiar function examples

  • Outlook

  • Q & A

Break

Second hour

  • Hands-on mock case studies

  • Q & A

Hello, I’m healthiar!

healthiar in numbers

  • 1 purpose - to help you quantify and monetize the attributable burden of disease

  • 3 dimensions - health, monetization, social aspects

  • 19 functions

  • > 1000 calculation pathways

  • 60000 (simple) assessments per second

healthiar overview

healthiar overview

healthiar overview

Example: attribute_health() with RR

Refresher - Burden of disease with relative risk

Refresher - Burden of disease with relative risk

Refresher - Burden of disease with relative risk

Refresher - Burden of disease with relative risk

attribute_health() with RR 1/3

Goal: attribute COPD cases to PM2.5 air pollution exposure

results_pm_copd <- attribute_health(
  erf_shape = "log_linear", # erf = exposure-response function
  rr_central = 1.369, 
  rr_increment = 10,  # μg / m^3
  exp_central = 8.85, # μg / m^3
  cutoff_central = 5, # μg / m^3
  bhd_central = 30747 # bhd = baseline health data (here: COPD incidence)
) 

Tip

healthiar comes with some example data that start with exdat_ that allow you to test functions.

results_pm_copd <- attribute_health(
  erf_shape = "log_linear",
  rr_central = exdat_pm$relative_risk, 
  rr_lower = exdat_pm$relative_risk_lower,
  rr_upper = exdat_pm$relative_risk_upper,
  rr_increment = 10, 
  exp_central = exdat_pm$mean_concentration,
  cutoff_central = exdat_pm$cut_off_value,
  bhd_central = exdat_pm$incidence
) 

attribute_health() with RR 2/3

attribute_health() outputs two lists (“folders”)

  • health_main contains the main results

  • health_detailed contains more detailed results and other assessment details in sub-lists (“sub-folders”)

Tip

Other healthiar functions also create ..._main & ..._detailed folders

results_pm_copd$health_main

Go to the Environment tab in RStudio …

… and click on a variable to “open” it.

Alternatively, you can use View(results_noise_ha), which has the same effect.

impact_rounded impact pop_fraction erf_ci rr exp bhd
3502 3501.962 0.1138961 central 1.369 8.85 30747

Example: compare() two scenarios

compare() two scenarios 1/2

Goal: compare attributable health impacts in two (policy) scenarios

  1. Use attribute_health() to calculate burden of scenarios 1 & 2
scenario_1 <- attribute_health(
  erf_shape = "log_linear",
  rr_central = 1.118,
  rr_increment = 10,  
  exp_central = 8.85, # EXPOSURE 1
  bhd_central = 25000
)
scenario_2 <- attribute_health(
  erf_shape = "log_linear",
  rr_central = 1.118,
  rr_increment = 10,
  exp_central = 6, # EXPOSURE 2
  bhd_central = 25000
)
  1. Use compare() to compare scenarios 1 & 2
results_comparison <- compare(
  output_attribute_scen_1 = scenario_1,
  output_attribute_scen_2 = scenario_2,
  approach_comparison = "delta" # or "pif" (population impact fraction)
)

compare() two scenarios 2/2

results_comparison$health_main
impact impact_rounded impact_scen_1 impact_scen_2 bhd exp_category exp_length exp_type exp_scen_1 exp_scen_2
731.5956 732 2349.958 1618.362 25000 1 1 population_weighted_mean 8.85 6

Example: attribute_health() across multiple geographic units

attribute_health() across multiple geographic units 1/3

Goal: attribute disease cases to PM2.5 exposure in multiple geographic units, such as municipalities, provinces, countries, …

results_iteration <- attribute_health(
  geo_id_micro = c("Zürich", "Basel", "Geneva", "Ticino", "Jura"), 
  geo_id_macro = c("German","German","French","Italian","French"),
  erf_shape = "log_linear",
  rr_central = 1.369,
  rr_increment = 10, 
  cutoff_central = 5,
  exp_central = c(11, 11, 10, 8, 7),
  bhd_central = c(4000, 2500, 3000, 1500, 500)
)

Here the we want to determine results by canton and then aggregate them by language region ("German", "French", "Italian")

attribute_health() across multiple geographic units 2/3

results_iteration$health_detailed$results_raw
geo_id_micro geo_id_macro impact_rounded
Zürich German 687
Basel German 429
Geneva French 436
Ticino Italian 135
Jura French 30

attribute_health() across multiple geographic units 3/3

results_iteration$health_main
geo_id_macro impact_rounded erf_ci exp_ci bhd_ci
German 1116 central central central
French 466 central central central
Italian 135 central central central
French 466 central central central

Tip

The main output contains aggregated results if available, or disaggregated results if no aggregation ID was provided

Example: summarize_uncertainty()

summarize_uncertainty() 1/4

If we add 95% confidence intervals to attribute_health()

results_pm_copd <- attribute_health(
  erf_shape = "log_linear",
  rr_central = 1.369, rr_lower = 1.114, rr_upper = 1.664,
  rr_increment = 10,
  exp_central = 8.85, exp_lower = 8, exp_upper = 10,
  cutoff_central = 5,
  bhd_central = 30747, bhd_lower = 28000, bhd_upper = 32000
) 

… then the folder health_detailed contains all different uncertainty combinations:

results_pm_copd$health_detailed$results_raw
erf_ci exp_ci bhd_ci impact_rounded rr exp bhd
central central central 3502 1.369 8.85 30747
lower central central 1252 1.114 8.85 30747
upper central central 5474 1.664 8.85 30747
central central lower 3189 1.369 8.85 28000
lower central lower 1140 1.114 8.85 28000
upper central lower 4985 1.664 8.85 28000
central central upper 3645 1.369 8.85 32000
lower central upper 1303 1.114 8.85 32000
upper central upper 5697 1.664 8.85 32000
central lower central 2765 1.369 8.00 30747
lower lower central 980 1.114 8.00 30747
upper lower central 4356 1.664 8.00 30747
central lower lower 2518 1.369 8.00 28000
lower lower lower 892 1.114 8.00 28000
upper lower lower 3967 1.664 8.00 28000
central lower upper 2877 1.369 8.00 32000
lower lower upper 1020 1.114 8.00 32000
upper lower upper 4533 1.664 8.00 32000
central upper central 4468 1.369 10.00 30747
lower upper central 1616 1.114 10.00 30747
upper upper central 6911 1.664 10.00 30747
central upper lower 4069 1.369 10.00 28000
lower upper lower 1471 1.114 10.00 28000
upper upper lower 6294 1.664 10.00 28000
central upper upper 4651 1.369 10.00 32000
lower upper upper 1682 1.114 10.00 32000
upper upper upper 7193 1.664 10.00 32000

summarize_uncertainty() 2/4

Goal: perform a Monte Carlo simulation to obtain summary uncertainty confidence intervals that combine uncertainty in different input data

Step 1: Use attribute_health() with uncertainty in ≥ 2 function arguments

results <- 
  attribute_health(
    erf_shape = "log_linear",
    rr_central = 1.369, 
    rr_lower = 1.114, 
    rr_upper = 1.664,
    rr_increment = 10,
    exp_central = 8.85, 
    exp_lower = 8, 
    exp_upper = 10,
    cutoff_central = 5,
    bhd_central = 30747, 
    bhd_lower = 28000, 
    bhd_upper = 32000) 

Step 2: summarize_uncertainty() to obtain the confidence intervals

results <- 
  summarize_uncertainty(
    output_attribute = results,
    n_sim = 100
  )

summarize_uncertainty() 3/4

results$uncertainty_main
impact_ci impact impact_rounded
central_estimate 3136.993 3137
lower_estimate 1537.385 1537
upper_estimate 5676.227 5676

summarize_uncertainty() 4/4

results$uncertainty_detailed$attribute_by_sim_disaggregated
sim_id impact rr_central exp_central bhd_central
1 5731.85 1.462496 10.16942 32127.59
2 3602.786 1.521621 7.911264 31318.83
3 4634.475 1.496932 8.774451 32812.33
4 2477.979 1.236879 9.054176 30007.59
5 3190.582 1.307062 9.020706 31257.61
6 3222.631 1.356192 8.512348 31754.01
7 3457.093 1.335157 9.418368 28834.7
8 2103.858 1.205937 8.744307 31070.1
9 2759.092 1.293466 8.644376 30822.27
10 2708.105 1.369631 8.044059 29659.24
11 4715.081 1.47613 9.283572 30688.8
12 2293.25 1.266201 8.123483 32268
13 5511.455 1.565692 9.435253 30564.39
14 1681.264 1.220164 7.876042 30226.44
15 1738.738 1.196227 8.373347 29645.66
16 2699.747 1.285776 8.672201 30618.67
17 3777.819 1.449824 8.667627 29662.75
18 3317.597 1.322765 9.046383 31000.83
19 2584.99 1.222993 9.448868 30176.27
20 3838.423 1.414973 8.716984 31710.98
21 6009.926 1.597893 9.638799 30756.55
22 5121.005 1.448811 9.864133 31034.71
23 2718.226 1.329628 8.382055 29591.6
24 2542.146 1.279985 8.520498 30542.15
25 2664.859 1.23451 9.037271 32682.35
26 1754.773 1.143348 9.277425 31509.66
27 4718.754 1.517077 8.941018 31152
28 5711.617 1.608763 9.226419 31374.23
29 4765.062 1.647204 8.661057 28534.14
30 3822.71 1.434593 8.834073 29583.24
31 2014.847 1.189432 8.774213 31791.97
32 5335.527 1.523692 9.280948 32342.64
33 3128.014 1.389939 8.1623 31633.04
34 2990.83 1.539019 7.437381 29982.26
35 3201.564 1.264264 9.618489 31191.94
36 2981.514 1.244678 9.641836 30861.85
37 1525.351 1.117177 9.518233 31236.91
38 1734.551 1.144898 9.332677 30461.48
39 5060.202 1.599173 8.62438 32339.76
40 5637.111 1.592109 9.465788 30058.56
41 4861.108 1.746363 8.229324 29502.68
42 2975.034 1.293848 8.888864 31207.69
43 2781.821 1.242423 9.565886 29482.27
44 1847.698 1.186965 8.53828 31400.12
45 2579.354 1.262943 8.848151 30021.75
46 4492.341 1.493217 8.894856 31072.59
47 3236.645 1.35796 8.623589 30839.8
48 2791.669 1.293426 8.720863 30578.25
49 2951.856 1.409966 8.048208 29688.17
50 5411.126 1.567719 9.337891 30536.99
51 3995.623 1.446675 8.702427 31268.44
52 3660.505 1.399679 8.826942 30316.47
53 2863.505 1.257943 9.117803 31757.71
54 2879.483 1.357981 8.300835 29972.25
55 4967.785 1.451873 10.0093 29159
56 3512.789 1.340506 9.137586 30763.13
57 2388.841 1.240475 8.876587 29806.86
58 4056.082 1.498765 8.591182 29989.67
59 3997.468 1.368704 9.493419 30390.03
60 4253.34 1.414422 9.236188 31137.12
61 2406.876 1.237385 8.963959 29726.87
62 2931.834 1.287537 8.873243 31440.49
63 3048.935 1.260693 9.421512 31316.64
64 857.4647 1.078443 8.817682 30172.27
65 4142.95 1.373768 9.598045 30495.46
66 2750.125 1.236879 9.452916 30447.9
67 3111.815 1.318392 8.854116 30793.45
68 2951.179 1.26288 9.158913 31903.06
69 4421.306 1.471209 8.885705 31737.16
70 3267.894 1.308333 9.261365 30199.23
71 3145.971 1.356827 8.746336 29122.17
72 3984.806 1.474773 8.590024 30608.95
73 4842.724 1.435846 9.57238 31765.57
74 4374.633 1.502303 8.985174 29217.65
75 4704.338 1.46052 9.738521 28631.75
76 2906.343 1.341438 8.341227 31089.41
77 4181.643 1.549701 8.636472 28396.42
78 2777.18 1.282869 8.75227 31122.64
79 2721.38 1.253935 9.026258 31250.85
80 3036.078 1.298591 9.138926 29620.38
81 4519.879 1.483334 9.085931 30376.04
82 2349.412 1.215117 9.081671 30732.38
83 2044.165 1.17841 9.042366 31836.65
84 3899.213 1.398355 8.905416 31769.13
85 4555.26 1.452661 9.383542 30169.96
86 3722.001 1.444415 8.48674 30931.49
87 2782.397 1.323915 8.536976 29449.77
88 5009.033 1.521841 9.204557 30948.65
89 1550.685 1.151462 8.629601 31075.18
90 3526.874 1.314995 9.476953 30568.3
91 2545.529 1.243599 8.98853 30565.73
92 2279.704 1.22658 8.776856 30709.3
93 3767.188 1.394531 9.16243 29141.71
94 3209.05 1.377354 8.523845 30078.43
95 2023.652 1.181774 8.956203 31649.47
96 2836.351 1.235014 9.633116 30443.77
97 961.7481 1.101541 8.349209 30175.87
98 4894.979 1.596689 8.666183 31050.85
99 2353.573 1.310017 7.882021 31433.46
100 3034.848 1.383293 8.134663 31381.57

Example: attribute_health() with RR and user-defined ERF

attribute_health() with RR and user-defined ERF 1/2

Goal: attribute COPD cases to air pollution exposure by applying a user-defined exposure response function (e.g. MR-BRT curves from Global Burden of Disease study)

  • Any function of the form intercept + a x c^1 + b x c^2 + …

  • Any other (non-linear) function of the function type, as obtained from e.g. splinefun() or approxfun()

results_pm_copd_mr_brt <- attribute_health(
  exp_central = 8.85,
  bhd_central = 30747,
  erf_eq_central = splinefun(
    x = c(0, 5, 10, 15, 20, 25, 30, 50, 70, 90, 110),
    y = c(1.00, 1.04, 1.08, 1.12, 1.16, 1.20, 1.23, 1.35, 1.45, 1.53, 1.60),
    method = "natural")
)

attribute_health() with RR and user-defined ERF 2/2

Example: attribute_health() with RR and exposure distribution

attribute_health() with RR and exposure distribution

Goal: attribute COPD cases to air pollution exposure categories

results_pm_copd <- attribute_health(
  erf_shape = "log_linear",
  rr_central = 1.369, 
  rr_increment = 10,
  exp_central = c(5, 6, 7, 8, 9, 10),
  prop_pop_exp = c(0.2, 0.1, 0.2, 0.3, 0.1, 0.1),
  cutoff_central = 5,
  bhd_central = 30747 
)
results_pm_copd$health_main
impact_rounded pop_fraction exp prop_pop_exp rr_at_exp
2177 0.0707956 5, 6, 7, 8, 9, 10 0.2, 0.1, 0.2, 0.3, 0.1, 0.1 1, 1.03190649221009, 1.06483100866533, 1.09880603094837, 1.13386507701522, 1.1700427342623

Example: attribute_health() with AR

Refresher - Burden of disease with absolute risk

Refresher - Burden of disease with absolute risk

attribute_health() with AR 1/2

Goal: attribute cases of high annoyance (HA) to noise exposure

Exposure category Exposure mean Population exposed
exp < 55 4’268’785
55 ≤ exp < 60 57.5 387’500
60 ≤ exp < 65 62.5 286’000
65 ≤ exp < 70 67.5 191’800
70 ≤ exp < 75 72.5 72’200
75 ≤ exp 77.5 7’700
results_noise_ha <- attribute_health(
  approach_risk = "absolute_risk",
  erf_eq_central = "78.9270-3.1162*c+0.0342*c^2",
  exp_central = c(57.5, 62.5, 67.5, 72.5, 77.5),
  pop_exp = c(387500, 286000, 191800, 72200, 7700),
)

attribute_health() with AR 2/2

results_noise_ha$health_main
approach_risk exp_category impact
absolute_risk 1, 2, 3, 4, 5 174231.8
results_noise_ha$health_detailed$results_raw
exp_category exp pop_exp impact
1 57.5 387500 49674.594
2 62.5 286000 50788.595
3 67.5 191800 46813.105
4 72.5 72200 23657.232
5 77.5 7700 3298.314

Other healthiar functions

Other healthiar functions

  • attribute_lifetable() for YLL
  • attribute_health() for YLD (using the dw & duration arguments)
  • get_daly() as the sum of YLL and YLD
  • monetize() monetizes (attributable) health impacts
  • cba() performs a cost-benefit analysis (CBA)
  • multiexpose() considers exposure to 2 exposures (e.g. two different air pollutants) at the same time to quantify the attributable health impacts
  • prepare_exposure() calculates (population-weighted) mean exposure by combining spatial exposure geographic units data
  • socialize() determines burden attributable to differences in a social indicator
  • attribute_mod() modifies an existing attribute_health() assessment
  • prepare_mdi() creates the BEST-COST MDI (Multidimensional Deprivation Index)

Outlook

Outlook

  • healthiar is publicly available on this GitHub page
    • healthiar might be expanded in the future (if funding available)
  • Python version coming out later in 2026
  • Package submitted to CRAN

Q & A

Happy coding! : )

Appendix

Addtional function examples

For additional function examples please see the function documentations and the package vignette.

The next slides will show you how to access them.

healthiar in RStudio 1/3

Figure: Landing page of the healthiar package in RStudio, where you find the package vignettes and function documentation.

healthiar in RStudio 2/3

Figure: Intro to healthiar vignette (= detailed guide to the package)

healthiar in RStudio 3/3

Figure: Help page of attribute_health()

Export results

write.csv(x = results_pm_copd$health_main, file = "exported_results/results_pm_copd.csv")
save(results_pm_copd, file = "exported_results/results_pm_copd.Rdata")
write.xlsx(x = results_pm_copd$health_main, file = "exported_results/results_pm_copd.xlsx")

Exported to .xlsx format

Visualize results

Visualization is out of scope of healthiar. You can visualize in